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Abstract 

An analysis method based on a deformation (as opposed to damage ) approach has been developed to model the 
strain rate dependent,, nonlinear deformation of woven ceramic matrix composites, such as the Reinforced Carbon 
Carbon (RCC) material used on the leading edges of the Space Shuttle. In the developed model, the cliffeiences in 
the tension and compression deformation behaviors have also been accounted for. State variable viscoplastic 
equations originally developed for metals have been modified to analyze the ceramic matrix composites. To account 
for the tension/compression asymmetry in the material, the effective stress and effective inelastic strain definitions 
have been modified. The equations have also been modified to account for the fact that in an orthotropic composite 
the in-plane shear response is independent of the stiffness in the normal directions. The developed, equations have 
been implemented, into LS-DYNA through the use of user defined subroutines (UMATs). Several sample qualitative 
calculations have been conducted, which demonstrate the ability of the model to qualitatively capture the featuies of 
the deformation response present in woven ceramic matrix composites. 


Introduction 


In the analysis of the impact behavior of foam on the leading edge of the Space Shuttle, 
accurate modeling of the deformation response of the Reinforced Carbon Carbon (RCC) leading 
edge material, a woven ceramic matrix composite with a plain weave fabric pattern, is critical. 
However, this material has several characteristics in the macroscopic deformation response 
which make material modeling efforts challenging. The stress-strain curve for this material is 
nonlinear, and the deformation response changes with strain rate. Furthermore, the uniaxial 
tension and compression responses are different, which implies that on the macroscopic scale the 

hydrostatic stresses have an effect on the material response. 

In ongoing work at NASA Glenn Research Center, efforts to analyze this material using LS- 
DYNA material model 58, MAT_LAMINATED_COMPOSITE_FABRIC, have been reasonably 
successful (Carney et al 2004). However, there are several features of model 58 which indicate 
that examining alternative methods of simulating the material response would be beneficial. 
First of all, in model 58 strain rate dependence of the material behavior is not explicitly 
accounted for within the model. However, in preliminary tests on the RCC material the data 
have indicated that this material may indeed have a nontrivial rate dependence. While 
preliminary efforts are underway to investigate adding strain rate dependence to model 58, these 
efforts are incomplete at this time. Furthermore, in model 58 all of the nonlinearity of the 
deformation response is due to damage mechanisms instead of deformation mechanisms 
(Matzenmiller et al 1995). Providing an alternative modeling capability in which the 
nonlinearity in the deformation response is assumed to be due to deformation mechanisms, 
allowing for a decoupling of the deformation and damage/failure response, could be more 

representative of the actual physical behavior. This is particularly true in this case since in the 
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theory on which model 58 is based (Matzenmiller et al 1995), the assumptions used to develop 
the damage/failure criteria applied in the equations were based on behaviors representative of 
polymer matrix composites. These assumptions may not necessarily be absolutely valid for the 
ceramic matrix composite material which is used for the Space Shuttle leading edges. 

At NASA Glenn Research Center, efforts have been underway for several years to model the 
strain rate dependent, nonlinear deformation of polymers with the effects of hydrostatic stresses 
included. In these efforts, the Bodner-Partom viscoplasticity model based on the state variable 
approach (Bodner 2002), originally developed for metals, has been modified in order to 
incorporate the hydrostatic stress effects known to be present in polymers (Goldberg et al 2003). 
These variations include modifying the effective stress and effective inelastic strain definitions, 
based on an application of the Drucker-Prager yield criterion (Khan and Huang 1995), to account 
for the effects of the hydrostatic stresses. 

In this paper, the constitutive equations with an associative flow rule described above have 
been applied and modified in order to model the nonlinear, strain rate dependent, deformation of 
plain weave woven ceramic matrix composites, including RCC, in which the differences in the 
tensile and compressive response are accounted for within the equations. Due to the plain weave 
fabric geometry of the RCC, for the purposes of this study the material is assumed to be isotropic 
in the normal directions, even though RCC is a composite and not a homogenous material. 
However, since RCC is indeed an ortho tropic material on the macroscopic level, the shear 
response is assumed to be uncoupled from and independent of the normal response. The full 
theoretical development of the constitutive equations will be described in this paper, along with 
the procedures to be used to determine the material constants. The procedures used to implement 
the equations within LS-DYNA through the use of user defined subroutines will also be 
described. Finally, several sample calculations will be presented which will demonstrate the 
ability of the model to qualitatively capture the features of the deformation response of plain 
weave woven ceramic matrix composite materials. 

Theoretical Development 

To analyze the nonlinear, strain rate dependent deformation of plain weave woven ceramic 
matrix composites, the Bodner-Partom viscoplastic state variable model (Bodner 2002), which 
was originally developed to analyze the viscoplastic deformation of metals above one-half of the 
melting temperature and has been successfully used in such instances, has been modified. In 
state variable models, a single unified strain variable is defined to represent all inelastic strains 
(Stouffer and Dame 1996). Furthermore, in the state variable approach there is no defined yield 
stress. Inelastic strains are assumed to be present at all values of stress, the inelastic strains are 
just assumed to be very small in the “elastic” range of deformation. State variables, which 
evolve with stress and inelastic strain, are defined to represent the average effects of the 
deformation mechanisms. 

An important point to note is that frequently, in state variable models, each of the state 
variables and material constants have a fairly explicit link to specific deformation mechanisms, 
such as dislocation pileups in metals (Stouffer and Dame 1996). Here, due to the complex nature 
of the material deformation mechanisms in woven ceramic matrix composites the linkage 
between the state variables and material constants to specific deformation mechanisms is not 
nearly as well defined. For this work, an internal stress state variable is defined to represent 
phenomelogically the resistance to inelastic deformation, a state variable is defined to represent 
the effects of hydrostatic stresses and a state variable is defined in order to scale the shear 
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stresses to ensure a consistent effective stress. The material constants are characterized in order 
to correlate with the experimentally observed deformation of the material, and are not tightly 
linked to specific deformation mechanisms. However, the equations developed here employ a 
fairly simple formulation with a minimum of material constants that are fairly simple to 
characterize, and capture the main features of the deformation response. 

For this study, temperature and moisture effects have been neglected, as only room 
temperature data are currently available. All of the nonlinearity and strain rate dependence 
observed in the material is assumed in this work to be due to inelastic deformation, where in 
reality the nonlinearity is most likely due to a mixture of deformation and damage. Furthermore, 
as mentioned previously in reality plain weave woven ceramic matrix composites like RCC are 
orthotropic materials, as is typical of composites. However, due to the plain weave architecture 
of the fibers, the material can be assumed to be isotropic in the normal directions, with 
independent shear responses. 

Associated Flow Equation and Evolution Equations 

In order to derive the associated flow equation for the components of the inelastic strain rate 
tensor for plain weave woven ceramic matrix composites with hydrostatic stress effects (different 
responses in tension and compression) included, a procedure similar to that employed to derive 
the Prandtl-Reuss equations in classical plasticity is utilized (Khan and Huang 1995) (Stouffer 
and Dame 1996). The formulation employed in this study is based on that used by Pan and co- 
workers (Li and Pan 1990) (Chang and Pan 1997) and Hsu, et al (Hsu et al. 1999). First, an 
inelastic potential function based on the Drucker-Prager yield criterion (Khan and Huang 1995) 
is defined as 


f=J7i + acr lk (1) 

where a kk is the sum of the normal stress components (assuming the usual rules of indicial 
notation (Stouffer and Dame 1996)) and is equal to three times the hydrostatic stress. The 
variable a is a state variable which controls the level of the hydrostatic stress effects. The term 
“aa kk ” incorporates the effects of hydrostatic stresses into the inelastic potential function. The 
variable J 2 * is a modified version of the second invariant of the deviatoric stress tensor. The 
modification is required due to the fact that for the orthotropic composites examined here the 
shear responses are independent of the normal terms in the stiffness matrix. Therefore, in order 
to ensure a consistent effective stress, the shear stresses must be multiplied by a scaling state 
variable |3. The modified version of the deviatoric second invariant is given as follows 

K =-([K -<T^) 2 +(cT 22 -C7 33 ) 2 + (tT33-cr il ) 2 ]+/ ff k2+C I 'l 2 3+ cr 23] (2) 

6 

where Gy are the components of the stress tensor. Note that for the current work all three shear 
stress components are scaled by the same value of p. In actuality, each shear stress component 
may require a unique scaling factor. However, for the current work the model is implemented 
using thin shell elements, which results in the transverse shear stresses being approximate at best. 
Furthermore, at the current time only in-plane shear data is available to be examined. As a 
result, only the correct scaling of the in-plane shear stresses is considered significant for this 
work. The effects of possible errors in the scaling of the transverse shear stresses are assumed to 
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not be significant. The components of the inelastic strain rate tensor, e ! t] , are assumed to be 

proportional to the derivative of the potential function with respect to the components of the 
stress tensor, ay, as follows 





(3) 


where X is a scalar rate variable. By taking the specified derivative of the potential function, the 
result becomes 


9 f s; 

2^/77 
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(4) 


where 5y is the Kronecker delta and Sy* are the modified components of the deviatoric stress 
tensor. For the normal stress components, Sy is equal to Sy, the components of the deviatoric 
stress tensor. However, to maintain consistency, for the shear stress components the terms are 
scaled as follows 

(v 


where (3 is the shear scaling state variable discussed earlier. 

The next step in defining the flow rule is to define the effective stress and the effective 
inelastic strain rate. The effective stress, a e , is defined as follows 

= -Jif = fT + V3 gkt h . . (6) 

To determine the value of the scalar X in terms of the effective inelastic strain rate, £ l e , the 
principal of the equivalence of the inelastic work rate is employed. With the aid of Equations 3 
and 4, the inelastic work rate, W 1 , can be expressed by the following 


W = <J- £ = (7 

’] y y 




U = crJ! 


(7) 


After applying the effective stress definition given in Equation 6 and simplifying, it can be 
shown that 


i = V3 e' e (8) 

By substituting Equation 8 and Equation 4 into Equation 3, a definition for the components of 
the inelastic strain rate tensor can be given as 
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= -fit 


2^ 


+ cxS :: 


(9) 


By utilizing Equation 7, it can be shown that the effective inelastic strain rate can be defined as 



( 10 ) 


where e[ is the effective deviatoric inelastic strain rate and £ ! m is the mean inelastic strain rate. 

The inelastic shear strain rates need to be scaled by a constant scale factor k due to the fact that 
for an orthotropic material the shear terms in the stiffness response are independent of the normal 
terms. When k equals one, which would happen for isotropic materials, the inelastic strain rate 
definition given here matches the effective inelastic strain rate definition given by Pan and co- 
workers (Li and Pan 1990) (Chang and Pan 1997). 

By following a similar format to the Bodner-Partom model (Bodner 2002) (Stouffer and 
Dame 1996), the following definition is specified 


s 
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= D„ exp 
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2 »] 
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( 11 ) 


where D 0 and n are material constants. D 0 represents the maximum inelastic strain rate, and n 
controls the rate dependence of the material. Z is an isotropic state variable which represents the 
resistance to inelastic deformation (internal stress). By substituting Equation 1 1 into Equation 9, 
the final form of the flow equation is determined to be 


et = 2D„ exp 


v 2 n 




l^ + aSii 


(12) 


Note that the elastic components of the strain rate must be added to the inelastic strain rate to 
obtain the total strain rate. 

The rate of evolution of the internal stress state variable Z, the hydrostatic stress effect state 
variable a and the shear scaling factor state variable (3 are defined by the equations 


^ CO 
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nT 

'Tsh 
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(13) 

a = q(aj -a)e’ e 

(14) 

/3 = q(A-P)K 

(15) 


where q is a material constant representing the “hardening” rate, and Zi, oq and (3j are material 
constants representing the maximum values of Z, a and (3 respectively. The initial values of Z, a 
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and (3 are defined by the material constants Z 0 , a 0 and (3 0 . The term e[ in Equations 13,14 and 

15 represents the effective deviatoric inelastic strain rate, which was defined in Equation 10. An 
important point to note is that in the original Bodner model (Bodner 2002), the inelastic work 
rate was used instead of the effective inelastic strain rate in the evolution equation for the internal 
stress state variable. However, for this work the inelastic strain rate was deemed easier to work 
with from both computational and characterization points of view, particularly in the 
incorporation of hydrostatic stress effects. Furthermore, the effective inelastic strain rate has 
been used in other state variable constitutive models (Stouffer and Dame 1996). Since 
hydrostatic stress effects were not considered in the original Bodner model (Bodner 2002), the 
evolution equation for a is new to this work. Furthermore, since these equations were originally 
developed for homogenous isotropic materials and not composites, the evolution equation for p 
is also new to this work. The state variables a and p are assumed to evolve in the same manner 
as the state variable Z. By using this assumption the value of q used in Equations 14 and 15 can 
be shown to be the same as the value of q used in Equation 13. 

Determination of Material Constants 

The material constants that need to be determined include D 0 , n, Z 0 , Zi, a D , oq, p o , Pi, k and q. 
These values, along with the elastic modulus, shear moduli and Poisson’s ratio are the input 
parameters required for the model. The procedure to be used is summarized here. More details 
on the procedure as it was used for isotropic polymers can be found in Goldberg, et al. (Goldberg 
et al 2003). The constants are determined using tensile and compressive stress-strain curves 
obtained from a set of constant strain rate tests, along with a shear stress-shear strain curve 
obtained for at least one strain rate. Each tension and compression test is conducted at a different 
total strain rate. 

The first step in the process is to obtain the values of oq and a G . To accomplish this task, 
Equation 6 is used in combination with stress-strain data from constant strain rate uniaxial 
tension and compression tests. The values of these constants are assumed to be rate independent, 
so the results from only one strain rate need to be used to find the needed parameters. In 
practical application of the methodology, the uniaxial tension and compression tests used do not 
have to be at the exact same effective strain rate. As long as the effective strain rates from the 
two tests are approximately equal, the values obtained have been found to be valid, particularly if 
quasi- static strain rates are examined. Two stress values from the tensile and compressive stress- 
strain curves are required in order to find these two constants. First of all, the “saturation stress”, 
the stress level at which the stress-strain curve flattens out and becomes horizontal, is required. 
At this stress level, the inelastic strain rate is assumed to be equal to the total applied strain rate. 
If the actual stress-strain curve terminates before the saturation stress is reached (as often 
happens in compression), the curve must be extrapolated to obtain an approximate value of the 
saturation stress. Secondly, the stress level at which the tension and compression curves 
becomes nonlinear (the “yield” stress), is required. One important point needs to be considered 
when determining the stress level designating the onset of nonlinearity. Frequently, in woven 
ceramic matrix composites, the tensile stiffness across nearly the entire stress-strain curve tends 
to be significantly softer than the compressive stiffness, most likely due to the fact that in tension 
significant matrix microcracks develop almost immediately upon loading, which weakens the 
material. As a result, in cases where this phenomena occurs the initial portion of the 
compression curve is assumed to represent the “elastic” regime of the composite material, and in 
tension the onset of nonlinearity in the stress-strain curve is assumed to occur at a very small, 
almost negligible stress value. 
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Given these definitions, Equation 6 can now be used to determine the values of oci and a 0 . 
The primary assumption used at this point (and assumed implicitly in Equation 6) is that the 
effective stress at saturation under uniaxial tensile loading at a particular strain rate is equal to 
the effective stress at saturation under compression loading at the same equivalent strain rate. 
Likewise, the effective stress at the point the stress-strain curve becomes nonlinear under tensile 
loading is equal to the effective stress at the point the stress- strain curve becomes nonlinear 
under compression loading. Therefore, assuming the value of a at saturation is equal to oci, and 
the value of a at the point the stress-strain curve becomes nonlinear is equal to a G , the following 
equations are obtained for the case of having data from uniaxial tension and compression tests 

cr„(l + V3ar 1 )=cT, c (l-A/3a 1 ) (16) 

tT Bl ,(l + V3aJ=<T„ fc (l-V3»J (17) 

where a st and a sc are the tensile and compressive stresses at saturation, respectively, and c n it and 
G„ic are the tensile and compressive stresses at the point where the respective stress-strain curves 
become nonlinear. Therefore, by substituting the tensile and compressive stresses at saturation 
into this equation, followed by the tensile and compressive stresses at the point where the 
respective stress-strain curves become nonlinear, the required constants can be determined. At 
this point, by applying Equation 6 for the tensile and shear curves, using the appropriate 
“saturation” and “yield” stresses in tension and shear, the values of p o and Pi can be determined 
from the following equations 


(l + Ty 

(18) 

, (l + 43 a 0 ) = «j3/3 0 T nl 

(19) 


where G st and t s are the tensile and shear stresses at saturation, respectively, and a n it and x„i are 
the tensile and shear stresses at the point where the respective stress-strain curves become 
nonlinear. 

The values of D 0 , n, and Zi are characterized as follows. The value of D 0 is currently 
assumed to be equal to a value of 10 4 times the maximum applied strain rate, which correlates 
with the maximum inelastic strain rate. To detennine the values of n and Zi, first Equation 12 is 
simplified to the case of uniaxial tensile loading, leading to the following expression 


e 1 =2£> 0 exp 


V".l 


(1 + V3cr)( 
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V3| 


+ (X 


cr 


( 20 ) 


where e 1 is the uniaxial tensile inelastic strain rate in a constant strain rate tensile test, cr is the 
uniaxial tensile stress, and the remainder of the terms are as defined earlier. The case of uniaxial 
tensile loading is used to determine these constants since for woven ceramic matrix composites 
the tensile stress-strain curves obtained experimentally are more likely to display a defined 
“saturation” stress than the compression curves, which as shown below is crucial for determining 
the material constants. Furthermore, Equation 20 can be rearranged as follows 
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The required constants are determined from a set of tensile stress-strain curves obtained from 
constant strain rate tests. Each curve in this set is obtained at a different constant strain rate. 
Data pairs of the total strain rate (equal to the inelastic strain rate at saturation) and saturation 
stress values from each curve are taken. For each strain rate, the data values are substituted into 
Equation 21, after the natural logarithm of both sides of the equation is taken, and represent a 
point on a master curve. The number of points in the master curve equals the number of strain 
rates at which the tensile tests were conducted. A least squares regression analysis is then 
performed on the master curve. As can be shown from Equation 21, the slope of the best-fit line 
is equal to -2n. The intercept of the best-fit line is equal to 2n(ln (Zi)). 

To determine the value of Z 0 , Equation 21 can be used again. In this case, the value of the 
tensile stress where the stress-strain curve becomes nonlinear for a particular constant strain rate 
tensile test is used, as well as the value of the approximate inelastic strain rate when the stress- 
strain curve becomes nonlinear. The point where the stress-strain curve becomes nonlinear is 
defined as the approximate point where the curve appreciably deviates from a linear 
extrapolation of the initial data. The strain rate used in the constant strain rate test divided by 
100 was found by trial and error to approximate the inelastic strain rate at this point reasonably 
well. Using this data, Equation 21 is solved for Z, which is assumed to be equal to the value of 
Z 0 . Compression data can also be used, particularly for the case when the tensile stress at which 
the stress-strain curve becomes nonlinear is approximated due to the differences in the initial 
tension and compression responses as discussed above. In this case, Equations 20 and 21 need to 
be reformulated based on Equation 12 given that compressive stresses are being considered as 
opposed to tensile stresses. Another point to note is that using the data from the lowest strain 
rate test available has been found to give adequate values of Z 0 . However, the calculations can 
be made using data from all the available strain rates, and an average taken if required to obtain 
the value of the constant. 

To determine the value for q for Equations 13, 14 and 15, first Equation 13 is integrated for 
the case of pure tensile loading (again since for woven ceramic matrix composites the tensile 
curves are more likely to display a distinct saturation stress), resulting in the following relation 

Z=Z 1 -(Z 1 -Zjexp(- <?£,') (22) 


where sj is the effective inelastic strain, computed using Equation 10. At saturation, the value of 
the internal stress Z is assumed to approach Zi, resulting in the exponential term approaching 
zero. Assuming that saturation occurs when the following condition is satisfied 


exp 


r -q £ l A 

1 + V3 a 
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(23) 
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the equation is solved for q, where e./ is the inelastic tensile strain at saturation. Also, in 
Equation 23 the effective inelastic strain at saturation is simplified to the case of uniaxial tensile 
loading through the use of Equations 9 and 10. If the inelastic strain at saturation is found to 
vary with strain rate, the parameter q is computed at each strain rate and regression techniques 
are utilized to determine an expression for the variation of q. If Equations 14 and 15 are 
integrated, an expression similar to Equation 23 is obtained. At saturation, the values of a and P 
are assumed to approach oq and Pi, and equations identical to Equation 23 are obtained, which 
would lead to the same value for q. Therefore, identical values of q are used in Equations 13, 14 
and 15. 

To find the value of k, Equation 10 is used along with results from a uniaxial tensile test and a 
pure shear test. Specifically, the effective inelastic strain at saturation for a tensile test at a 
particular total strain rate (usually the lowest) is set equal to the effective inelastic strain at 
saturation from a pure shear test at the same approximate total strain rate, resulting in the 
following expression which can be solved for k 


1 + a / 3 a 



(24) 


where z} is once again the inelastic tensile strain at saturation and yy is the inelastic engineering 
shear strain at saturation. 


Finite Element Implementation 

The constitutive equations described above have been implemented into LS-DYNA through 
the use of the user defined subroutine (UMAT) option. The routines are designed to work with 
shell elements, therefore the equations are specialized for the case of plane stress in the normal 
directions. A composite is a thin structure, which justifies the use of shell elements. The use of 
shell elements allows greater flexibility in terms of defining a laminate of many composite plies, 
and in specifying the angle of orientation of the material. Furthermore, the use of shell elements 
provides greater flexibility and computational efficiency in terms of analyzing a complex 
structure. All of the material constants are entered into the input deck. The history variables for 
the finite element analysis are defined to include the deviatoric stresses, the hydrostatic stress, 
the values of the state variables and the inelastic strain rate tensor. To integrate the constitutive 
equations, both a standard Euler integrator (Kreyszig 1992) and a fourth order Runge-Kutta 
integrator (Kreyszig 1992) have been implemented within separate UMAT routines. The Euler 
routine is very simple to define and implement. However, in order to ensure accurate results, 
each time step must be broken up into ten mini- steps within the UMAT routine. The Runge- 
Kutta routine, while somewhat more complex to implement, does not require the time step 
division within the routine. Studies are still being carried out as to which integrator yields the 
most accurate answers with the greatest computational efficiency. 

In the UMAT routine, the material constants and the current values of the history variables 
are read into the routine. The values of the elastic stiffness tensor are computed based on the 
elastic modulus and Poisson’s ratio. Within the integration routines, for each integration step the 
effective stress is computed employing the current stress values and the current value of the state 
variable a using Equation 6. Next, the tenns in the inelastic strain rate tensor are determined 
using Equation 12. The effective inelastic strain rate is determined using Equation 10, and then 
the state variable rates are determined using Equations 13, 14 and 15. The appropriate 
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integration formula and time step increment is then used to compute updated values of the 
inelastic strain increments and the state variables. From the given total strain increments and the 
computed inelastic strain increments, along with the elastic constants, the values which are 
computed include the strain increments for the current step in the integration process, the stress 
increments and revised total stresses. In this process, the out of plane normal strain increments 
are calculated in a manner such that the plane stress condition is enforced, and the out of plane 
normal stress is set equal to zero. As is common for shell elements, a correction factor of 0.833 
is used in computing the transverse shear stresses. Revised values of the deviatoric and 
hydrostatic stresses determined from the computed stresses, and the integration process is 
continued as required. Finally, the history variables are written out, and the routine returns. 

Sample Calculations 

To demonstrate the ability of the model described above to qualitatively capture the features 
of the deformation response of woven ceramic matrix composites such as RCC, several sample 
analyses have been carried out using the current model. The tension, compression and shear 
stress-strain curves for a representative woven ceramic matrix composite with a plain weave 
fiber architecture were determined and are shown qualitatively in Figure 1. Specific numerical 
values for the stress-strain curve, the material constants that were utilized and comparisons to 
actual test data cannot be presented at this time due to formal NASA restrictions on the use and 
discussion of actual data on the RCC material. However, as can be seen from the presented 
stress-strain curves, the nonlinearity of the deformation response most likely present in actual 
experimental data is captured well. Furthermore, the presented computations demonstrate the 
ability of the constitutive equations to capture qualitatively the significant differences in the 
tension and compression responses using a single set of material constants. In Figure 2, tensile 
stress-strain curves at two different strain rates are generated. As can be seen in the figure, the 
model is able qualitatively to capture strain rate dependencies in the material response, which 
again is a factor in the actual deformation response of these materials. Future efforts will involve 
further exercising the model and comparing the computed results to publicly available test data. 


Conclusions 

The Bodner-Partom viscoplastic state variable model (Bodner 2002) has been modified in 
order to analyze the strain rate dependent, nonlinear deformation of woven ceramic matrix 
composites such as the Reinforced Carbon Carbon found in Space Shuttle leading edges. The 
differences in the tensile and compressive deformation responses have been accounted for 
through the modification of the effective stress and effective inelastic strain definitions given in 
the original equations. The procedures used to compute the shear response have also been 
modified to account for the fact that in an orthotropic composite the shear response is 
independent of the normal response. The material constants can be obtained based on tension, 
compression and shear tests. The constitutive equations have been implemented within LS- 
DYNA through the use of the user defined subroutine (UMAT) option. Sample calculations 
have been conducted which demonstrate the ability of the model to qualitatively capture the 
nonlinearity, strain rate dependence and tension/compression asymmetry observed in the 
deformation response of woven ceramic matrix composites such as Reinforced Carbon Carbon. 
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Figure 1: Computed tension, compression and shear deformation of a representative woven 

ceramic matrix composite material. 
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Figure 2: Effect of strain rate on the tensile deformation of a representative woven ceramic 

matrix composite material. 
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